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FOREWORD 


This is a progress report on the research project, "Analysis and 
Computation of Internal Flow Field in a Scramjet Engine," for the period 
ending December 31, 1986. The work was supported by the NASA Langley 
Research Center (Computational Methods Branch of the High-Speed 
Aerodynamics Division) through research grant NA6-1-423. The grant was 
monitored by Dr. A. Kumar and Mr. J.P Drummond of the High-Speed 
Aerodynamics Division. 



INTERACTION OF TRANSIENT RADIATION 
IN NONGRAY GASEOUS SYSTEMS 


S. N. Tiwari* and D. J. Singh'*’ 

Old Dominion University, Norfolk, VA 23508 

ABSTRACT 

A general formulation is presented to investigate the transient radiative 
interaction in nongray absorbing-emitting species between two parallel plates. 
Depending on the desired sophistication and accuracy, any nongray absorption 
model from the line-by-line models to the wide-band model correlations can be 
employed in the formulation to investigate the radiative interaction. Special 
attention is directed to investigate the radiative interaction in a system 
initially at a uniform reference temperature and suddenly the temperature of 
the bottom plate is reduced to a lower but constant temperature. The interac- 
tion is considered for the case of radiative equilibrium as well as for 
combined radiation and conduction. General as well as limiting forms of the 
governing equations are presented and solutions are obtained numerically by 
employing the method of variation of parameters. Specific results are 
obtained for CO, CO 2 , H 2 O, and OH. The information on species H 2 O and OH is 
of special interest for the proposed scramjet engine application. The results 
demonstrate the relative ability of different species for radiative 
interactions. 
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1. INTRODUCTION 

The field of radiative energy transfer in gaseous systems is getting an 
ever increasing attention recently because of its applications in the areas of 
the earth's radiation budget studies and climate modeling, fire and combustion 
research, entry and reentry phenomena, hypersonic propulsion and defense- 
oriented research. In most studies involving combined mass, momentum, and 

energy transfer, however, the radiative transfer formulation has been coupled 
mainly with the steady processes^"^^ and the interaction of radiation in 

transient processes has received very little attention. Yet, the transient 

approach appears to be the logical way of formulating a problem in a general 
sense for elegant numerical and computational solutions. The steady-state 
solutions can be obtained as limiting solutions for large times. 

A few studies available on radiative ineractions reveal that the 
transient behavior of a physical system can be influenced significantly in the 
presence of radiation^^”^^. Lick investigated the transient energy transfer 
by radiation and conduction through a semi-finite medium^^. A kernal 
substitution technique was used to obtain analytic solutions and display the 
main features and parameters of the problem. Steady and transient heat 
transfer in conducting and radiating planar and cylindrical mediums were 

analyzed in Refs. 13 and 14 according to the differential formulation. The 
analyses based essentially on the gray formulation provide some qualitative 
insight into the effect of absorption and emission on the transient 
temperature distribution in the gas. Ooornink and Hering^^ studied the 
transient radiative transfer in a stationary plane layer of a nonconducting 
medium bounded by black walls. A rectangular Milne-Eddington type relation 
was used to describe the frequency dependence of the absorption coefficient. 
It was found that the cooling of the layer initially at a uniform temperature 
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is strongly dependent on the absorption coefficient model employed. Larson 
and Viskanta^® investigated the problem of transient combined laminar free 
convection and radiation in a rectangular enclosure. It was demonstrated that 
the radiation dominates the heat transfer in the enclosure and alters the 
convective flow patterns significantly. The transient heat exchange between a 
radiating plate and a high- temperature gas flow was investigated by Melnikov 
and Sukhovich^^. Only the radiative interaction from the plate was 
considered; the gas was treated as a non-participating medium. It was proved 
that the surface temperature is a function of time and of longitudinal 
coordinate. Some other works on transient radiation and related areas are 
available in Refs. 18-22. 

The goal of this research is to include the nongray radiative formulation 
in the general unsteady governing equations and provide the step-by-step 
analysis and solution procedure for several realistic problems. The specific 
objective of the present study is to investigate the interaction of nongray 
radiation in transient transfer processes in a general sense. Attention, 
however, will be directed first to a simple problem of the transient radiative 
exchange between two parallel plates. In subsequent studies, the present 
analysis and numerical techniques will be extended to include the flow of 
homogeneous, nonhomogeneous, and chemically reacting species in one- and 
multi-dimensional systems. 
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2. BASIC THEORETICAL FORMULATION 


The physical model considered for the present study is the transient 
energy transfer by radiation in absorbing-emitting gases bounded by two 
parallel gray plates (Fig. 2.1). In general, Tj and T 2 can be a function of 
time and position and there may exist an initial temperature distribution in 
the gas. It is assumed that the radiative energy transfer in the axial 
direction is negligible in comparison to that in the normal direction. 

For radiation participating medium, the equations expressing conservation 
of mass and momentum remain unaltered, while the conservation of energy, in 
general, is expressed as^ 

P Cp ^ = div (k grad T) + pT ^ + p 4 ) - div (2.1) 


where p is dynamic viscosity, p is the coefficient of thermal expansion of the 
fluid and <t» is the Rayleigh dissipation function. For a semi-infinite medium 
capable of transferring energy only by radiation and conduction, Eq. (1) 
reduces to 


P c 


P 


dT _ _ dq 
■5T “ ” "Sy 


( 2 . 2 ) 


where q is the sum of the conductive heat flux ~ ~ (9T/9y) and the 


radiative flux q^. For the physical model where radiation is the only mode of 
energy transfer, the energy equation can be written as 


dT 


aqR 

dy 


(2.3) 


Use of this simplified equation is made to investigate the transient behavior 
of a radiation participating medium. 


For many engineering and astrophysical applications, the radiative 
transfer equations are formulated for one-dimensional planar systems. For 
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diffuse nonreflecting boundaries and in the absence of scattering, the 
expression for the total radiative flux is given, for a n-band gaseous system, 

t3yl>8,23 


qn(y) = e, - ep + I E / 7 

^ ^ i=l Aw, 0 ^“i “i ^ ""i 


■ ^ '"w. 7 <c,.U-y)]<lz}dw (2.4) 

y 1 1 1 


where 


‘ - ®2co. 

11 111 1 

Equation (2.4) is in proper form for obtaining the nongray solutions of 

molecular species. In fact, this is an ideal equation for the line-by-line 

and narrow-band model formulations. However, in order to be able to use the 

wide band models and correlations, Eq. (2.4) is transformed in terms of the 

correlation quantities 

q„ (!) • ej - *01 "o1 ^7 %((!-!')]«' 


1 - 1 

- / ^2^.(0 A( C|u^. (C‘-0]d^*} (2.5) 

C 1 

where 

I = u/Uq = y/L; V = u'/u^ = z/L; A = A/A^; 
u « (S/A_,) Py; . (S/A^) PL; PS = 

It should be noted that and F 2 (^, (2.5) represent the values at 

the center of the ith band and A'(u) denotes the derivative of ^(u) with 
respect to u. Upon performing the integration by parts, Eq. (2.5) can be 
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expressed in an alternate form as 

n £ 


23 


q^(C) = e^ - 62 + A^. {/ Cde^ (5')/dC'] A. [| )]dC' 


+ / Cde^ (C')/d?'3 A. [j Uq. (?*-5)]dC'} 


( 2 . 6 ) 


A direct differentiation of Eq. (2.6) provides the expression for the 
divergence of radiative flux as 

dqn(^) . n E 

- 3 T- = 7 ,^/of “01 ' 

- 1 

a: [f u^^.(E-5')]dE’ - Cde^_(5')/dE’] x 

a;. [| Uoi(E'-C)]dE'} (2.7) 

Equations (2.5) through (2.7) are the most convenient equations to use when 
employing the band-model correlations in radiative transfer analyses. 

By defining (|) (E.t) = T(E,t)/T^ with Tq representing some constant 
reference temperature, Eqs. (2.2) and (2.7) can be combined to yield the 
energy equation (for the general case of simultaneous conduction and 
radiation) in nondimensional form as 

(E»t)/5t = R w - •» S {/ (p . (E* »t) X 
^ ^ 1=1 0 

- 1 1 

A( [| u^. (E-C')]dE‘ - (C',t) X 

a; [| u^^. (E'-?)]dE'} 


( 2 . 8 ) 
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where 

R = k t /(p C L^) 
m p 

(k^.(5,t) = {PS.(T)[5 e^.U,t)/SC]/(p Cp yt,^)} 

The time t in Eq. (2.8) is defined as t* = t/t^ with t^ representing some 

characteristic time scale of the physical problem; however, for the sake of 

convenience, the asterisk is left out here as well as in further developments. 

From the definitions of 4i(C,t) and .(C,t), it should be noted that Eq. 

0)1 

(2.8) is a nonlinear equation in T(C,t). Equation (2.8), therefore, 

represents a general case of the transient energy transfer by radiation and 
conduction between two semi-infinite parallel plates. The nondimensional 

parameter R in Eq. (2.8) is analogous to the Fourier number. For R = 0, Eq. 

(2.8) reduces to the case of pure transient radiative energy transfer. The 
initial and boundary conditions for Eq. (2.8) will depend on the conditions of 
the specific physical problem. 


2.1 Specific Application 

As a special case, it is assumed that the entire system initially is at 
the fixed (reference) temperature Tq. For all time, the temperature of the 

upper plate is maintained at the constant temperature equal to the reference 
temperature, i.e., T£ = Tq. The temperature of the lower plate is suddenly 
decreased to a lower but constant value, i.e., T^ < T£. The problem, 
therefore, is to investigate the transient cooling rate of the gas for a step 
change in temperature of the lower plate. 

In many radiative transfer analyses, it is often convenient (although not 
essential) to employ the relation for the linearized radiation as 


1 


1 


(d 


1 w 


(T 




(2.9) 
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where again the subscript i refers to the ith band such that is the wave 
number location of the band and represents the temperature of the reference 
wall which could be either Tj or T 2 . For the special case considered, since 
we are interested in investigating the transient behavior of the gas because 
of a step change in temperature of the lower plate, is taken to be equal to 
Tj. It should be pointed out that for a single-band gas, the linearization is 
not required because the temperature distribution can be obtained from Eq. 
(2.8) and the radiative heat flux can be calculated from Eqs. (2.5) and (2.6). 
However, for the case of multiband gases and for systems involving mixtures of 
gases, it is convenient to employ the linearization procedure in order to use 
the information on band model correlations. The following definitions are 
useful in expressing the governing equations in linearized forms: 

9 = (T-Tj)/(T 2 -Tj) (2.10a) 

"li ' '’'V*' s’ '^11 • 

n 

Nj = (P tjQ Cp) , Kj = Kj. (2.10c) 

”li = (y- P «^p) »li ‘ «li = 

Ml = (yi P Cp) Hj. (2.10e) 


“^li ‘^oi ' ^li ' 


^ oi "li = PL 'li 


(2.10f) 


where Hj, K^, Nj and Mj represent the values of H, K, N and M evaluated at the 
temperature Tj. As explained in Refs. 1 and 8, these quantities represent the 
properties of the gaseous medium. 
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The initial and boundary conditions for the physical problem considered 
are 

9(C.O) = 1; 9(0, t) = 0; 9(1, t) = 1 (2.11) 

It is important to note that the boundary conditions given in Eq. (2.11) are 
applicable only to the case of simultaneous conduction and radiation energy 
transfer. The cases of transient radiation and radiation with conduction are 
treated separately in the following sections. 

2.2 Transient Radiation 

By employing the definitions of Eqs. (2.9) and (2.10), Eq. (2.8) is 
transformed to obtain a convenient form of the energy equation for the 
transient radiation case as^^ 

_ S9(c,t) . 3 ^ ^ ae(^'.t) 

5t 2 . . ^li bV 

1=1 0 

T.r3 r.n^r. ^ ae(5'.t) ^ 

J • a;' ” 

a;c| (2.12) 

The parameters in Eq. (2.12) are and Uq. For a given gas, the parameters 
are the gas pressure and the temperature of the lower wall. As pointed out 

earlier, the boundary conditions given in Eq. (2.11) are not applicable to Eq. 

(2.12) because this equation does not require a boundary condition. Thus, in 
this case, the temperature of the medium adjacent to a surface differs from 
the surface temperature. This is because the temperature of the medium 

adjacent to a surface is affected not only by the surface but also by all 

other volume elements and surfaces. The radiation slip method is a means of 
accounting for such temperature jumps and this is discussed in Ref. 1. 
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For the case of transient radiative interaction, the nondimensional 
radiative heat flux is defined by 


Q(^,t) = qp(5,t)/[ej(t) - e2(t)] 


(2.13) 


By employing the definitions of Eqs. (2.9), (2.10), and (2.13), relations for 
the radiative flux, as given by Eqs. (2.5) and (2.6), are expressed as 

n 5 

f 

0 


Q(^,t) = 1- (3/8 a Tp u^^. {/ 0(^',t) x 


A! cl u^.(C-E')]dC +/ Cl - e(C'.t)] X 

1 c U I ^ 

a;.c| u^j.(C'-E)]dC‘> 


(2.14a) 


and 


Q(C,t) = 1 - (1/4 o Tj^) t H .{j^ 

^ i=l 0 


a;| u^.(C-C')]dC* - / 


1 


bQ(V ,t) 

— sr~ 


A.c| u^i(C'-C)]dC'} (2.14b) 

It should be pointed out that Eq. (2.14a) is a convenient form for the 
optically thin and general solutions while Eq. (2.14b) is useful for solutions 
in the large path length limit. Once the solutions for 9(^,t) are known 
from the energy equation, the appropriate relations for the heat flux can be 
obtained from Eqs. (2.14). It should be noted that the quantity Hj/(aT^) in 
Eqs. (2.14) is nondimensional. 
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2.3 Radiation and Conduction 

For this case where conduction heat transfer takes place simultaneously 
with radiative transfer, the energy equation is given by Eq. (2.8). Thus, a 
combination of Eqs. (2.8)-(2.10) results in 

ae _ „ s ^0 _ 3 " „ ae(c',t) „ 


A. - / 


dO(V ,t) 
— 5T“ 


A] [| u^.(C* - C)]d^'} (2.15) 

Since the presence of conduction implies continuity of temperatures at the 
boundaries, the boundary conditions for Eq. (2.15) are those given in Eq. 
(2.11). The quantity Nj^/R can be expressed as N = N^/R = (PL^/k)K 2 . The 
nondimensional parameter N denotes the relative importance of radiation versus 
conduction in the gas. For particular values of P and L, it is actually the 
dimensional gas property N 2 /(RPL^) = Kj/k that represents the relative 
importance of radiation versus conduction. 

For the case of no radiation, Eq. (2.15) becomes 

II .rJ (2.16) 

The separation of variables results in a general solution of Eq. (2.16) as 

e = exp(-X^ Rt) (Bj sin cos \ 5) (2.17) 

p 

where \ is the separation parameter. The particular solution of Eq. (2.17) 
can be obtained by satisfying the boundary conditions. Alternately, by 
defining a similarity variable n = ^//4Rt, Eq. (2.16) can be written as 


12 


^.2.^=0 (2.18) 

The solution of Eq. (2.18), with initial and boundary conditions given by Eq. 
(2.11), is found to be 


0(5,t) = (2//i) erf Ti (2.19) 

This solution is applicable for relatively large separations between the 
plates. 

In the case of simultaneous conduction and radiation heat transfer, the 
nondimensional heat flux is defined as 

Q = Cqj.(C,t) + q^(C,t) ]/[ej^(t) - e2(t)] (2.20) 

Alternately, this can be expressed as 

Q = C(d0/a^) + Q (2.21) 

where 

C = k(Tj - T2)/{L[ej(t) - e2(t)]} 

The expression for Q in Eq. (2.21) is obtained either from Eq. (2.14a) or 
from Eq. (2.14b). It should be noted that the relation for 0 needed in Eq. 
(2.21) comes from the solution of Eq. (2.15). 
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3. METHOD OF SOLUTIONS 

The solution procedures for the transient radiation and radiation with 
conduction cases are presented in this section. The procedure for the case of 
transient radiative energy transfer is given first and this is followed by the 
procedure for combined radiation and conduction. 


3.1 Transient Radiative Equilibrium 

For the general case of transient radiative equilibrium, the temperature 
distribution is obtained from the solution of the energy equation, Eq. (2.12). 
Once 9(i^,t) is known, the radiative heat flux is calculated by using the 
appropriate form of Eq. (2.14). Before discussing the solution procedure for 
the general case, however, it is desirable to obtain the limiting forms of 
Eqs. (2.12) and (2.14) in the optically thin and large path length limits and 
investigate the solutions of resulting equations. 


3.1.1 Optically Thin Limit 

In the optically thin limit A(u) = u and A'(u) = 1, and therefore, Eq. 
(2.12) reduces to^»*^»^^ 


+ 3 Nj 9(5, t) - I Nj = 0 


(3.1a) 


From an examination of Eq. (3.1a) along with the definitions given in Eq. 
(2.10), it is evidient that in the optically thin limit the temperature 
distribution in the medium is independent of the 5- coordinate for the case 
of pure radiative exchange. This is a charateristic of the optically thin 
radiation in the absence of other modes of energy transfer. Thus, Eq. (3.1a) 
can be wri tten as 


d9(t) 

dt 


+ 3 N^ 


9(t) 


0 ; 0 ( 5 . 0 ) = 1 


(3.1b) 


14 


Since gas properties are evaluated at known reference conditions, is 

essentially constant, and solution of Eq. {3.1b) is found to be 

9(t) = |[1 + exp(- 3 Njt)] (3.2) 

In the optically thin limit both forms of Eq. (2.14) yield the same final 
relation for the radiative flux as^^ 

Q(C,t) = 1 - [3/(8 otJ)] (PLKj) C(l-C) + 

+ (2C - 1) e(t)] (3.3) 

3 

It should be pointed out that in Eq. (3.3) the quantity (PLKj/aTj^) is 

nondimensional . The relation for 9{t) in Eq. (3.3) is obtained from Eq. 
(3.2). Thus, evaluation of the temperature distribution and radiative heat 
flux in the optically thin limit does not require numerical solutions. 

3.1.2 Large Path Length Limit 

In the large path length limit (i.e., for u^^- >> 1 for each band), one 
has A(u) = = Jln(u), A'(u) = 1/u, and A''(u) = - 1/u . Thus, in this limit, 
Eq. (2.12) reduces to^*®*^^ 

99(5, t) _ _ „ 99(5', t) d5' M /i^ 

An analytical solution of Eq. (3.4) may be possible, but numerical solution 
can be obtained quite easily. 

In the large path length limit, Eqs. (2.14a) and (2.14b) reduce 


respectively to 
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Q(5,t) = 1 - (l/cr4T^) H^[/ e(C'.t) 

0 


and 


- / 


dC' 


TFT7 


T - r ] 


(3.5a) 


Q(C,t) = 1 - (l/4aT^) Z H {/ X 

1 . , n 0^ 

1=1 0 


^n[| Uo.(C-C')]dC' +/^M|^x 


^n[| u^.(C'-C)]dC'} (3.5b) 

The expressions for dimensionless radiative heat flux from or to the wall are 
obtained by setting C = 0 in Eqs. (3.5). 

3.1.3 Numerical Solutions of Governing Equations 

The general solutions of Eqs. (2.12) and (3.4) are obtained numerically 
by employing the method of variation of parameters. For this, a polynomial 
form for 9(5, t) is assumed in powers of 5 with time dependent coefficients as 

n 

0(5, t) = E cjt) 5 (3.6) 

m=o 

By considering only the quadratic solution in 5, and satisfying the boundary 
conditions of Eq. (2.11), one finds 

9(5, t) = 5^ + g(t) (5-5^) (3.7) 

where g(t) represents the time dependent coefficient. At t = 0, a combination 
of Eqs. (2.11) and (3.7) yields the result 


g(0) = (l-5^)/(5-5^) 


(3.8) 
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Equations (3.7) and (3.8) are used to obtain specific solutions of Eqs. (2.12) 
and (3.4). Once the temperature distribution is known, the expressions for 

heat flux are obtained by numerically integrating the corresponding equations. 
The entire numerical procedure is described in detail in Ref. 23. The 
numerical procedure is similar for higher order solutions in C of Eq. (3.6), 
but computational resources required are considerably higher. 


3.2 Radiation with Conduction 


Equations (2.15), (2.20) and (2.21) are the appropriate equations to 
describe the combined process of conduction and radiation energy transfer. 
The energy equation, Eq. (2.15), reduces in the optically thin limit to 


&t 


+ 3N^ 9 


3 

I 


2 

d e 


(3.9) 


A closed form solution of this equation may be possible but a numerical 
solution can be obtained easily. In the large path length limit, Eq. (2.15) 

reduces to 


dt 


= R 




f 99 dC' 

Jq w tftt 


(3.10) 


The solution of this equation is obtained numerically. 

The general solution of Eq. (2.15) is obtained by using the method of 
variation of parameters as discussed in Sec. 3.1 and Ref. 23. For a quadratic 
solution in C, a combination of Eqs. (2.15) and (3.7) results in 

g‘ (t) + J^(^) g(t) = J 2 (^) 


(3.11a) 
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where 

o n C 

J.(C) = 2 R C(C) +4 C(C) 2 N.. {/ (1-2I') a: Cb.(^-C')]dC' 


1 

- / (1-2C') a;. [b.(?'-5)]d5'} (3.11b) 


n I 

JM = 2 R C(C) - 3 C(C) Z N,. {/ V a;, [b . (5-C' ) ]dC‘ 
^ 1=1 ^^0 ’ ^ 


1 

- / V A’ [b.(?’-C)]dC'} (3.11c) 


^1 " ^ (3. lid) 

and n represents the number of vibration-rotation bands. The solution of Eq. 
(3.11a) is found to be^^ 

g(t) = [g(0) - J 2 (C)/Ji(?)l exp[-Jj(?)t] + J2(5)/Ji(C) (3.12) 

where g(0) is given by Eq. (3.8). The integrals in and J 2 are evaluated 
numerically by following the procedure discussed in Ref. 23. With g(t) known 
from Eq. (3.12), Eq. (3.7) provides the solution for the temperature distribu- 
tion. Once 9(C,t) is known, the heat transfer is calculated by using either 
Eq. (2.20) or Eq. (2.21). 

In the optically thin limit, the expression for and J 2 » as defined by 
Eqs. (3.11b) and (3.11c), reduce to 

0^(5) = 2 R C(5) + 3 Nj (3.13a) 

J2(C) = Oil) [2 R - - 1/2)N^] 


(3.13b) 
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In the large path length limit, the expressions for and are obtained as 


J (O = C(V {2 R / C(1-2C')/(C-C')M'} 

^ 0 


(3.14a) 


1 

Jp(C) = C(C) (2 R - 3 M, / Cc7(C-C')]dC’} (3.14b) 

^ 0 


Equations (3.13) and (3.14) are useful in obtaining numerical solutions for 
the limiting cases. 


19 


4. PHYSICAL CONDITIONS AND DATA SOURCE 

For the physical problem considered (Fig. 2.1) four specific absorbing- 
emitting species were selected for an extensive study; these are CO, CO 2 , OH 
and H 2 O. The species CO was selected because It contains only one fundamental 
vibration-rotation (VR) band and all spectral Information are easily available 
In the literature. It Is a very convenient gas to test the numerical 
procedure without requiring excessive computational resources. Species OH and 
H 2 O are the primary radiation participating species for the pressure and 
temperature range anticipated In the combustor of the scramjet engine. 
Species CO 2 , and combinations of CO 2 and H 2 O are Important absorbing-emitting 
species In many other combustion processes. 

In radiative transfer analyses. It Is essential to employ a suitable 
model to represent the absorption-emission characteristics of specific species 
under Investigation. Several llne-by-llne, narrow-band, and wide-band models 
are available to model the absorption of a VR band (Refs. 7-11). However, It 
Is often desirable to use a simple correlation to represent the total 
absorption of a wide band. Several such correlations are available In the 
literature (Refs. 7-11). The relative merits of these correlations are 
discussed In Ref. 11. In this study, the correlation proposed by Tien and 
Lowder^ Is employed and this Is given by 

A(u) = In (uf(p) 'I + 1> (4.1) 

where 

f(p) = 2.94[1 - exp(- 2.60 p)] 


and p represents the line structure parameter. 
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The spectral information and correlation quantities needed for CO, CO 2 

and H 2 O were obtained from Refs. 7-9. The spectral data for OH were obtained 
« 

from Refs. 24 and 25, and correlation quantities are developed in Appendix A. 
The specific VR bands considered for each species are: CO {4.7ti funda- 

mental), OH (2.8|i fundamental), CO 2 (15|a, 4.3h and 2.7^), and H 2 O (20ii, 
6.3|i, 2.7n, 1.87 |a, and 1.38fi). The radiative properties of important 

species are provided in Appendix A. 

It is important to consider the variation of thermodynamic and transport 
properties with temperature and pressure. This information is available for 
important species in Refs. 26-28 and essential information for the present 
study is provided in Appendix B. 

For the specific problem considered, the dependent variables are 0, Q, 
and Q, and independent variables are C and t. The parameters for general- 
solutions are Tj, T 2 /TJ, Uq, and t^. For the radiative equilibrium case, 0 
and Q depend only on t and in the optically thin limit and on t, and 

in the large path length limit. For the case of combined radiation and 

conduction, 0 and Q depend on t, and N = N^/R in the optically thin limit 

and on C, t, and M = M^/R in the large path length limit. Information on 

radiative ability of various species in the optically thin and large path 
length limits is available in Ref. 8 and Appendix C. The parameters for 
specific solutions for different species are Tj^, T 2 /T 1 , P, L, and tj„. 

Extensive results, therefore, can be obtained by varying these parameters. 
For parametric studies, however, only certain values of pressure, temperature, 
and plate spacing were selected and results were obtained for the general as 
well as limiting cases. Unless stated otherwise, specific results were 
obtained for T 2 = 2 Tj^. 
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5. RESULTS AND DISCUSSION 

Results have been obtained for different radiation participating species 
for both cases, the radiative equilibrium and radiation with conduction. The 
computer program used for numerical solutions is provided in Appendix D. It 

should be realized at the outset that, according to the physics of the 
problem, the gas initially is at a high temperature Tj = T 2 . At t = 0, the 

temperature Tj is lowered to a constant value. The energy exchange then 

occurs and the gas cools down in time until a steady-state conditon is 
reached. At this time, a certain temperature profile is established and a 
fixed amount of energy exchange occurs irrespective of the time. The rate of 
cooling of the gas layer, therefore, is dependent on the nature of the 
participating species and on the physical parameters of the problem. 

Some limiting solutions that are independent of any participating species 
are presented first in Figs. 5.1 and 5.2 for the radiative equilibrium case. 
In Fig. 5.1, the temperature distribution in the channel is plotted as a 
function of the optically thin parameter Nj^ for different times. Similar 

results are illustrated in Fig. 5.2 (with solid lines) for even higher values 
of Nj. These results show an exponential decay with time reaching the steady- 
state value of 9=1/2 for t->®. The temperature distribution for the large 
path length limit is shown in Fig. 5.2 (with broken lines) as a function of 

the large path length parameter and for different times. Although the 
numerical values are entirely different, these results also show the 

exponential decay with time and reach the limiting value of 9 = 1/2 for 
t ->■ ®. It should be noted that while the optically thin solutions are 
independent of the ^-coordinate, the large path length solutions do depend on 
5 and they have been obtained for C = 0.5. In the case of simultaneous 
radiation and conduction, both optically thin and large path length solutions 
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Fig. 5.1 Temperature variation with the optically 
thin parameter N, . 
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Limiting solutions for the radiative 
equilibrium case. 
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for temperature distribution depend on These results, however, can be 

expressed in terms of the radiation-conduction parameter N=Nj/R in the 
optically thin limit and M = Mj/R in the large path length limit. 

The radiative equilibrium temperature distribution for CO are shown in 
Fig. 5.3 for three different characteristic times and for P = 1 atm, L * 10 
cm, and T^ = 500 K. For small t^, t* = t/t^ becomes large and, therefore, 9 
varies slowly with t*; the reverse is true for large t^ values. This trend is 
evident clearly from the simplest case of the optically thin solution given by 
Eq. (3.2). Similar trends in results were observed also for other species for 
different values of P, T^, and L. Thus, to demonstrate typical transient 
trends, other results presented in this study were obtained for an 
intermediate value of the characteristic time t„, = 0.00001 sec. 

The centerline temperature variations with time are illustrated in Fig. 
5.4a for CO and in Fig. 5.4b for OH. General and limiting solutions are shown 
for pure conduction, radiation, and radiation with conduction for P = 1 atm, 
Ty, = 500 K, and L = 10 cm. It is noted that for both gases the optically thin 
solutions approach the steady-state conditions faster than the large path 
length and general solutions. For the physical conditions considered, the 
energy is transferred faster by conduction than by radiation, and the steady- 
state conditions are reached earlier by the combined radiation and conduction 
process. Although OH is a relatively better heat conducting gas, CO is seen 
to be more effective in the radiative transfer. For the same physical 
conditions as in Figs. 5.4, the radiation and radiation with conduction 
results are compared for CO, OH, H 2 O, and CO 2 in Fig. 5.5. It is seen that H 2 O 
is most effective and OH is least effective in transferring the radiative 
energy. The ability of a gas to transfer radiative energy depends on the 
molecular structure of the gas, band intensities and physical conditions of 



Fig. 5.3 Temperature variation for transient radiation for 
different characteristic times (CO, P = 1 atm, 

T = 500 K, L = 10 cm and C = 0.5). 
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Fig. 5.4a Centerline temperature variation with time for CO 
with P = 1 atm, T, , = 500 K, and L = 10 cm. 









with P = 1 atm, T , = 500 K, and L = 10 cm. 



Fig. 5.5 Comparison of centerline temperature variation with time for 
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the problem. Thus, H 2 O with five strong VR bands is a highly radiation 
participating species and the steady-state conditions are reached quickly for 
H 2 O than for other species. However, CO with one fundamental band is seen to 
be a better radiating gas than CO 2 with three VR bands. This is because for 
the given physical conditions, the optical thickness of CO 2 is sufficiently 
large and in the large path length limit CO 2 is relatively less effective in 
transferring the radiative energy (Ref. 8 and Appendix C). Further results 
for CO and OH are illustrated in Figs. 5.6 for different wall temperatures. 
It is seen that while radiation is less effective than conduction at T^ = 500 
K, it is highly effective at T^ = 1,000 K. This, however, would be expected 
because radiation becomes considerably important at higher temperatures. The 
steady-state condition is reached quicker for T^ = 1,000 K than for T^ = 500 
K. In fact for the characteristic time considered (tjj, = 0.00001 sec.), the 
steady-state condition is reached quickly for all species for temperatures 
higher than T^^ = 1,000 K. Results for the pure radiation case are illustrated 
in Figs. 5.7 for CO and OH for L = 1 cm and 10 cm. It is seen that while the 
general and large path length solutions depend on the plate spacing the 
optically thin solutions are independent of the spacing. This fact was 
pointed out earlier in the method of solution. In the presence of other modes 
of energy transfer, the optically thin solutions also depend on the plate 
spacings. As would be expected, for the same physical conditions, the steady- 
state condition is reached quicker for the lower plate spacing. 

The temperature variations within the plates are shown in Figs. 5.8-5.11 
for different species and for P = 1 atm and L = 5 cm. In the absence of 
molecular conduction, temperature jumps (radiation slips) occur at the 
boundaries and, therefore, the general solutions for the case of radiative 
equilibrium are not presented. However, general as well as limiting solutions 


Radiation and conduction 






Fig. 5.6b Centerline temperature variation with time for OH with P 







Fig. b.8a Temperature variation for combined radiation and conduction 
for CO with P = 1 atm, T . = 1000 K, and L = 5 cm. 
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Fig. 5.8b Temperature variation for CO with P = 1 atm, T =1000 K, and L = 5 cm. 
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Fig. 5.9a Temperature variation for combined radiation and conduction 
for OH with P = 1 atm, T . = 500 K, and L = 5 cm. 



Fig. 5.9b Temperature variation in optically thin limit for OH 
with P = 1 atm, T, . = 500 K, and L = 5 cm. 








Fig. 5.10a Temperature variation for combined radiation and conduction 
for HpO with P = 1 atm, T =500 K, and L = 5 cm. 
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Fig. 5.11 Comparison of temperature variation for combined radiation 
and conduction for P = 1 atm, T. . = 500 K, and L = 5 cm. 
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are presented for the case of radiation with conduction. It is noted, in 
general, that for the case of radiation with conduction, the steady-state 
conditions are reached for all species at t > 0.1 and for > 500 K. For 
the case of pure radiation, the steady-state conditions are reached at 
relatively longer times. The optically thin results are seen to be independent 
of the ^-coordinate for the case of radiative equilibrium and are seen to 
vary slowly in the central portion of the plates for the case of radiation 
with conduction. This is because, in this limit, the gas interacts directly 
with the boundaries and conduction is predominant near the walls. Specific 
results for CO are illustrated in Figs. 5.8 for T^^ = 1,000 K. For the case of 
radiation with conduction, general and limiting solutions are compared in Fig. 
5.8a; and for both cases, the radiative equilibrium and radiation with 
conduction, limiting solutions are compared in Fig. 5.8b. The steady-state 
results for pure conduction are also shown in Fig. 5.8b for comparative 

purposes. The results demonstrate the typical trends for limiting and general 
solutions, i.e., a lower temperature gradient implies a higher rate of energy 
transfer. Specific results for OH are illustrated in Figs. 5.9 for T^^ = 500 
K. General and limiting solutions are shown in Fig. 5.9a; and limiting 
solutions are compared in Figs. 5.9b and 5.9c. For the case of radiation with 
conduction, the limiting and general solutions are seen to be in good 
agreement for all times (Fig. 5.9a). This is because for the conditions of 
the illustrated results, conduction dominates the energy transfer process in 
OH. The typical trends in results for the optically thin and large path 
length limits are shown in Figs. 5.9b and 5.9c, respectively. Figure 5.9b 

clearly shows that for all times the radiative equilibrium results are 
independent of the ^-coordinate in the optically thin limit. Figure 5.9 

shows that at earlier times the rate of energy transfer is higher in the 
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presence of conduction. Specific results for H 2 O are illustrated in Figs. 

5.10 for the case of radiation with conduction. It is seen clearly that the 
rate of cooling is signficantly higher in the large path length limit (Fig. 
5.10a), and the steady-state conditions are reached at relatively longer times 
for lower values (Fig. 5.10b). For the case of combined radiation and 
conduction, a comparison of results for different species is shown in Fig. 

5.11 for t = 0.01 and 0.1. The results for t = 0.1 essentially correspond to 
the steady-state conditions. For t = 0.01, the variation in temperature is 
seen to be relatively small between 5=0.2 and 0.9. The centerline 
temperature is found to be the lowest for H 2 O, and this is followed by OH, CO, 
and CO 2 . However, it is noted that OH is very effective in transferring the 
net energy in comparison to the other species. As discussed earlier, this is 
mainly due to relatively higher conductive ability of OH at T^ = 500 K. 

The centerline temperature distributions are shown in Figs. 5.12-5.15 for 
different gases as a function of the spacing between the plates. In most 
figures, results are presented for both cases, the radiative equilibrium and 
radiation with conduction. In selected figures, results for the case of pure 
conduction are included also for comparative purposes. For a particular gas, 
specific results are presented for various times to demonstrate the radiative 
nature of the gas under different pressure and temperature conditions. 

The results for CO are presented in Figs. 5.12 for different cases. For 
P = 1 and T^ = 500 K, the results illustrated in Fig. 5.12a show that the time 
required to reach the steady-state condition increases with increasing plate 
spacings. For a particular plate spacing, the centerline temperature is lower 
for the case of radiation with conduction than for pure radiation for all 
times. For P = 1 atm and T^^ = 1,000 K, results presented in Fig. 5.12b show 
that the large path length solutions are closer to the general solutions for L 
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> 20 cm; and the results for pure radiation and radiation with conduction are 
identical for t > 0.5. The centerline temperature variations are shown in 

Fig. 5.12c for t = 0.5 * 500 K, and different pressures. It is noted that 

while the heat transfer by conduction is insensitive to the change in 
pressure, the radiative heat transfer is strongly dependent on it. The rate 
of radiative interaction increases with increasing pressure until the large 
path length limit is reached for sufficiently large values of L. For the case 
of pure radiation, the results for P = 0.1 atm differ considerably from other 
results. This is due to use of the Tien and Lowder's correlation which is 
suitable only at relatively higher pressures (Ref. 11). The centerline 

temperature variations are shown in Fig. 5.12d for t = 0.5, P = 1 atm, and 
different values of T^^. As would be expected, both conductive and radiative 
interactions increase with increasing temperatures, although the increase in 
radiative transfer is comparatively higher. It should be noted that for T^^ = 
300 K, T 2 = 2 T^,, = 600 K, for T^^ = 500 K, T 2 = 1,000 K, and so on. Thus, for 
a higher value of T,^ = T^, the energy interactions occur at a sufficiently 
large temperature difference between the upper and lower plates. At these 
temperatures, if the plate spacing is small, the energy is transferred quickly 
and the steady-state condition is reached at relatively shorter times. This 
fact was pointed out also in the discussion of results of Figs. 5.7. 

The centerline temperature variations for OH are illustrated in Figs. 
5.13 for different conditions. The results presented in Fig. 5.13a for P = 1 
atm and T^ = 500 K show the similar trend as CO in Fig. 5.12a, although the 
extent of energy transfer by simultaneous radiation and conduction is 
relatively higher. This is because at T^ = 500 the energy transfer in OH is 
dominated by the conduction heat transfer. General and limiting solutions for 
radiative equilibrium are shown in Fig. 5.13b and for radiation with 
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conduction in Fig. 5.13c. These results clearly demonstrate the typical 
radiative interaction trends for different times. The results show that the 
optically thin solutions are independent of the plate spacing in the case of 
pure radiation but depend on the spacing when molecular conduction is 
included. The large path length results are seen to be valid only for large 
values of L for the case of pure radiation (Fig. 5.13b), but they appear to be 
valid in the entire range for the case of radiation with conduction (Fig. 
5.13c). The results for pure conduction, pure radiation, and conduction with 
radiation are illustrated in Fig. 13d for P = 1 atm and T^^ = 1,000 K. For 

this temperature, the results for pure radiation and radiation wih conduction 

are found to be identical. This indicates that at higher temperatures, OH 
becomes a highly radiation participating gas. The results for variation of 
the centerline temperature for OH with pressure and temperature are given in 
Figs. 5.13e and 5.13f and they show the same general trend as the results for 
CO shown in Figs. 5.12c and 5.12d. 

Extensive results of 0(5 = 0.5) versus L have been obtained for H 2 O and 
CO 2 for different conditions and some of these are illustrated in Figs. 5.14. 
In general, these results show similar trends as exhibited by the results for 
CO and OH but the extent of radiative interactions is entirely different. For 
example, a comparison of results presented in Fig. 5.14a for H 2 O with the 
results of Fig. 5.12a for CO and Fig. 5.13a for OH for identical conditions 

reveals that H 2 O is a highly radiation participating gas even for shorter 
times. For H 2 O, results presented in Fig. 5.14b demonstrate that the large 
path length solutions are closer to the general solutions for both cases, the 
radiative equilibrium and radiation with conduction. The centerline 
temperature variation with L is shown for H 2 O in Figs. 5.14c and 5.14d for 

different pressure and temperature. As epxected, radiative effects are seen 
to be higher at higher pressure and temperature. 



Centerline temperature variation with L for CO 
(P = 1 atm, and T, = 500 K). 
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Fig. 5.12c Centerline temperature variation with L for CO ( T = 500 K 
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Centerline temperature variation with L for OH 
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Fig. 5.13b Centerline temperature variation with L for pure radiation 
(OH, P = 1 atm, and T = 500 K). 







Centerline temperature variation with L for combined radiation 
and conduction (OH, P = 1 atm, and T = 500 K). 
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Fig. 5.13d variation with L for OH (P = 1 atm. 







Centerline temperature variation with L for OH at different 
pressures and T = 500 K. 



Fig. 5.13f Centerline temperature variation with L for OH at different 
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The centerline temperature variations are compared for different gases in 
Fig. 5.15 for P = 1 atm, = 500 K, and t = 0.05. For the case of radiative 

equilibrium, it is noted that OH is the least effective and H 2 O is the most 

effective gas in transferring the radiative energy for plate spacings greater 
than two centimeters. When molecular conduction is included, OH becomes more 
effective because of its relatively higher conductive ability. These points 
were noted also in earlier discussions. The story, however, can be entirely 
different for other physical conditions because of the radiative/conductive 
nature of participating species (Ref. 8 and Appendix C). This fact is 

partially evident from the steady-state results, for the case of combined 
radiation and conduction, presented in Fig. 5.16 for two different 

temperatures, T^ = 300 K and 500 K. For example, for T^^^ = 300 K and L = 10 
cm, the temperature values for CO and CO 2 are about the same, for H 2 O it is 
lower, and for OH it is the lowest; however, for plate spacing greater than L 
= 20 cm, the trend is entirely different. Also, it should be noted that the 
steady-state (t = 0.5) results for Ty^ = 500 K in Fig. 5.16 show different 
trend than the results for the same temperature in Fig. 5.15 for t = 0.05. 
Thus, in order to predict the relative ability of a gas for radiative 
interactions, it is very important to consider the exact physical conditions 
for all species. These predictions may not be applicable if physical 

conditions of the problem are changed. 

Extensive results for the variation in heat transfer can be presented 
analogous to the variation of temperature for different conditions. However, 
this should not be necessary because the heat transfer variation follows the 
trend of the temperature variation in a reverse manner. If the temperature 
differences are higher, the rate of heat transfer will be higher and the 
steady-state conditions will be reached at earlier times. The results for 
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Fig 5.15 Comparison of centerline temperature variation with L for 
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heat transfer variations have been obtained for selected conditions and these 
are discussed here briefly. 

For P = 1 atm, the results for Q and Q are illustrated, as a function 
of t*, in Figs. 5.17-5.20 for different species. The results for t* = 0.0 - 
1.0 are shown in Fig. 5.17a for H2O. However, it is found that for the t^ 
value selected in this study, the steady-state conditions are reached in most 
cases at about t* = 0.2. Consequently, other results (Figs. 5.17b, 5.17c, 
5.18-5.20) are presented only in the range of t* = 0.0 - 0.2 to demonstrate 
the rate of cooling at different temperatures. As would be expected, the 
results show that for a given plate spacing the gas layer reaches the steady- 
state condition faster at higher values of because of stronger radiative 
interactions. It should be noted that the rate of energy transfer increases 
with time for a gas layer closer to the upper wall (C = 0.75) and decreases 
with time for a gas layer closer to the lower wall (C = 0.25) until the 
steady-state conditions are reached. The rate of cooling is entirely 
different if the plate spacing is changed (Fig. 5.17c). From a comparison of 
results of Figs. 5.17-5.20, it is noted that the trend and rate of energy 
transfer are different for different species. This, however, would be 

expected because of the relative ability of different species to participate 
in the radiation-conduction interaction process (Appendix C). 


RAD & COND Gas: HjO 

RADIATION tm = 0. 00001 sec 
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Fig. 5.17c Variation in heat flux with time for H 20 (p = 1 atm, L = 10 cm 
and t* = 0.0 = 0.2). 







RAD & COND Gas: OH 

RADIATION tfn = 0. 00001 sec 

^ = 0.25 P = 1atm 

1^=0. 75 L = 5 cm 
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Fig. 5.20 Variation in heat flux with time for OH (P = 1 atm, 
and t* = 0.0 - 0.2). 
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6 . CONCLUSIONS 

The problem of transient radiative interaction in nongray absorbing- 
emitting species has been formulated in a general sense such that sophisti- 
cated absorption models can be used to obtain accurate results if desired. 
Results have been obtained for the special case of radiative interactions in a 
plane gas layer bounded by two parallel plates when the temperature of the 

bottom plate is suddenly reduced to a lower but constant temperature. The 
energy transfer by pure radiation, and by simultaneous radiation and 
conduction were considered and specific results have been obtained for CO, 
CO2. H2O, and OH by employing the Tien an Lowder's correlation for band 
absorption. It is noted that the extent of radiative interaction is dependent 
on the nature of the participating species and parameters Tj, T2/T2, P, L, and 
tjj,. The steady-state conditions are reached at relatively longer times for 
radiative equilibrium than for radiation with conduction. For a particular 
value of P and T^, the time required to reach the steady-state condition 

increases with increasing plate spacing. For a fixed plate spacing, the 

energy is transferred quickly for higher Tj values because of large 
temperature differences between the plates. The rate of radiative interaction 
increases with increasing pressure until the large path length limit is 

reached. The radiative equilibrium solutions are found to be independent of 
the plate spacing in the optically thin limit. In the case of simultaneous 
radiation and conduction, both optically thin and large path length solutions 
depend on the y-location between the plates. At moderate temperatures, OH is 
a poor radiating but better heat conducting gas. For most conditions, H2O is 
found to be highly radiation participating species, and the steady-state 
conditions are reached quickly for H2O than for other species. The 

information on radiative interactions for OH and H2O is useful in the analysis 
of the flow field in the scramjet engine. 
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APPENDIX A 


SPECTRAL INFORMATION AND CORRELATION QUANTITIES 
FOR IMPORTANT INFRARED BANDS 


Spectral information and correlation quantities for important infrared 
bands are available in [7-9]. The exponential band model correlation 
quantities for important bands of CO, CO 2 , H^O, and CH^ are given in Table 
A.l. For some species, revised data are now available in [93- 

In developing correlation quantities for different bands, it is 
important to have information on the integrated band intensity S(T^) and the 
band width parameter A^(T). The band intensities of many molecular bands 
are available in the literature [7-9, 24, 25] but this is not the case witn 
the band width parameter. Different relations for A^(T) are proposed by 
Edwards et al. in references cited in [9]. From a critical evaluation of 
different relations available in the literature for Aq(T) and after personal 
communications with Dr. D.K. Edwards, it was decided to adopt the following 
two relations: 


Ao(T) 


1 r2 ^ 4kTBj 1/2 

2 4"^ ^ he ^ 


(A.l) 


A,(T) 


0.9 


^ ^kTB^ 1/2 
4 he 


(A.2) 


where B is the rotational constant of lower level, c is the speed of light, 
h is the Planck's constant and k is the Boltzmann's constant. The value of 
B is different for different molecules. Equations (A.l) and (A.2) are 


A-1 


essentially the same except Eq. (A.l) will result in a coefficient of 0.707 
instead of 0.9. It is suggested to use Eq. (A. 2) in determining approximate 
relations for A^(T) for all molecules whose values are not available in the 
literature. Information on the rotational constants for different 
molecules is available in the literature (for example, see Ref. 24). It is 
suggested to use the equivalent value for the rotational constant of 
polyatomic molecules. 

It should be noted that in Eqs. (A.l) and (A. 2), B, c, h, and k are 
constants and do not depend on the temperature. Thus, Eq. (A. 2) may be 
expressed as: 


Ag(T) = CONST (T)^^^ (A. 3a) 

where 

CONST = 0.9 (i) (illfi) (A. 3b) 

4 he 

and represents the equivalent rotational constant. By evaluating Eq. 

(A. 3a) at a reference temperature The value of Aq(T^^^) can be 

determined and, therefore, Eq. (A. 3a) may be expressed alternately as: 

Equation (A. 4) is a convenient form to compare its results with experimental 
values. 

2 3 2 

By noting that r (_) = [r (3/4)] and substituting values for c, h, 

4 

and k, Eq. (A. 3) can be expressed as: 

Aq(T) = 1.59313 (BgT)^/2 


(A.5) 
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where and Bg have units of cm ^ and T is in degrees Kelvin. For a 
particular gas, A^ obtained from Eq. (A. 5) and then Eq. (A. 4) 

can be used to determine A^(T) at other temperatures. For example, for CO 
the rotational constant is 1.931 cm”^ and at a reference temperature of 
300 K, A^CT^gf * 300 K) = 38.344 cm This compares very well with the 
experimental value of 38.1 given if Refs. 7 and 9 and presented in Table 
A.l. Similarly, for the 4.3 y band of CO 2 , the equivalent rotational 
constant is 0.3906 cm'^ and, therefore, A^(T^g^ = 300 K) = 17.246 cm”^; the 
experimental value of 19.9 cm”^ given in Table A.l for this band is slightly 
higher. 

Spectral Information for OH 

For the fundamental band of OH, the following information is obtained 
from Ref. 25: 


Band center, = 3570 cm’^, 

-2 -1 

Band strength at STP, S(T^) = 110 cm atm 

Also, from Ref. 24 the information on equilibrium rotational constant for OH 
is obtained as: 

A 17.355 cm‘^ 

X 18.871 cm"^ 
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Thus, it is suggested to use a value for the equilibrium rotational constant 
for OH as = 18 cm Using this, a value of (T^gf= 300 K) is found 

as 

Aq (T^g^: = 300) = 1.59313 (BgT^g^)^^^ = 117.0707. (A.6) 

Considering the value of = 117 cm \ the relation for A^(T) for OH 

is given as: 

Aq(T) = 117 (T/300)^^^. (A. 7) 

By knowing S (T^), and A^(T), other required spectral information for 
OH can be evaluated. 


Table A.l Exponential Band Model Correlation Quantities** 












A -6 


1. Correlation quantities are based on the results of reference [9], The 
intensity of the band marked with * was taken from the reference [7]. 

2. Notations: K^(T) = (T/300)^^^ K 2 (T) = (300/T)^'^^, 

= [^1 ^(T)/K^(T)] X 10'^ . «2 = *^3 ^T)/[^2(T)Ki(T)] 

h = 6.625 X 10 erg-sec, C = 2.998 x 10^^ cm/sec, 
k = 1.380 X 10 erg/K , hC/k = 1.44 cm - K 



<)) 2 {T) = {1 - exp [( - hC/KT ) 

{[1 - exp ( - hCu^/kT)] [1 - exp (-hCoj 3 /kT) ]}‘^ , ^> 3 ( 1 ) = 1 + 
3/2 


0.053 (T/lOO) 



6 . 


For H^O, ( 0 ^ = 3652 cm"^, = 1595 cm"^, = 3756 cm 

'**''1 ''2 ''3 = {1 - exp [ - hC(v^ co^ + o) 3 )/kT]} x 

{[1-exp ( -hCu^/kT)] [l-exp(-hC( 02 /kT] [l-exp( -hCu^/kT]} 

(T) = exp [ - 17.6 (T/100)*^^^] 



APPENDIX B 


THERMODYNAMIC AND TRANSPORT PROPERTIES 
OF SELECTED SPECIES 

It is important to consider the variation of thermodynamic and trans- 
port properties of various species with temperature and pressure. Quite 
often, the variation in properties with the pressure is not as crucial as 
with the temperature. The information on variation of different properties 
is provided here. 

The information on variation of the thermal conductivity with temper- 
ature is obtained from Ref. 27 and this is expressed as: 

where 

T = Temperature, K 

= thermal conductivity at T^ = 273 K, k cal/m-hr-K 
n = constant as given in Table B1 

In order to be consistent with the units used in the present work, Eq. (B.l) 
is expressed as: 


X 


4 

10 X, 


(36) (23.889) 


(T/273)" 


erg/sec-cm-K 


(B.2) 
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B-2 


For Eq. (B.2) the value of is obtained from Table B.l; the units for 
are shown in the table. For example, thermal conductivities of CO and CO, 
are expressed as: 

CO: X = (11.627865289 x 200) (T/273)°‘® 

=2325.570579* ( (TW/273.0)**0.8) 

CO^ : X = (11.627865289 x 128) (T/273)^*^^ 

= 1488.365171* ( (TW/273.0)**1.23) 

Thus, for Eq. (B.2), the tabulated values for X^ should be used without 

4 

dividing by the factor 10 . 

Thermal conductivities of other species can be calculated in a similar 
manner. For species not listed in the Table B.l, values should be obtained 
from Ref. 27 or other sources such as Ref. 26. For higher tmperatures, 
values available in Ref. 26 should be used. Some of the values used in the 
present study are listed in Table B.2. 

The relations for the constant-pressure specific heat for different 
ideal gases are available in Ref. 28 and these are given in Table B.3 for 


(B.3) 


(B.4) 


CO, OH, CO^, and H^O. 













Table B.2 Thermal Conductivity of Selected Species 

erg/(cm-sec-K). 


Temp., K 


1,000 


2,000 


Molecule 

OH 

CO 

1 

CO 2 

H 2 O 

(4879.71) 

2507.82 

(2674.22) 

1671.43 

(1820.47) 

1738.04 

(2925.31) 

(6993.13) 

3773.77 

(3938.09) 

3133.02 

(3339.63) 

3703.67 

(4980.15) 

(11504.56) 

6570. 51 

7349.02 

10325. 76 


(6888.51) 

(6716.93) 

(11588.26) 

(20276.33) 

11439. 93 

17238.37 

28803. 58 


(11730. 56) 

(11822.63) 

(26302. 73) 

1 


* Values in parenthesis are from Ref. 26 
















APPENOIX C 


RADIATIVE ABILITY OF SELECTED SPECIES IN OPTICALLY THIN 
AND LARGE PATH LENGTH LIMITS 

The quantities, M^ and N^ defined in Eqs. (2.10) represent the 
interaction parameters for radiation versus conduction heat transfer. As 
pointed out in Secs. 2.3 and 4, the dimensional gas property Nj^/k represents 
the relative importance of radiation versus conduction in the optically thin 
limit. Similarly, the gas property Mj/(RL) = Hj^/k represents the relative 
importance of radiation versus conduction in the large path length limit. 

The radiative properties and Hj^ are defined, respectively, by Eqs. 

(2.10c) and (2.10e). 

The quantities and K^/k were calculated for a number of gases and 
the results are illustrated in Figs. C.l. The results presented in Fig. 

C.la show the radiative ability of different species in the optically thin 
limit. The results presented in Fig. C.lb demonstrate the relative 
importance of radiation versus conduction in the optically thin limit. 
Similarly, the results for Hj and Hj^/k were calculated to investigate the 
interaction in the large path length limit and these are illustrated in 
Figs. C.2a and C.2b. A comparison of results presented in Figs. C.l and C.2 
shows a considerable difference in the radiation-conduction interaction for 
the optically thin limit as opposed to the large path length limit. For 
example, in the optically thin limit CO 2 possesses a large radiation 
interaction relative to the other gases, while the reverse is true in the 
large path length limit. On the other hand, just the opposite trend is 
observed for H 2 O. 


C-1 




Fig. C.la Radiative ability of different species in the optically thin 









in the large 





Hi/k = M/L, cm" 
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Temperature , K 


Fig. C.2b Interaction parameter for 
large path length limit. 



APPENDIX D 


COMPUTER PROGRAM FOR TRANSIENT RADIATIVE TRANSFER IN NONGRAY GASES 

Following the numerical procedure of Sec. 3, a general computer program 
is developed to calculate the transient radiative transfer in nongray gases 
between two parallel plates having constant wall temperature. The program 
can be used to calculate the conduction and radiation heat transfer 
separately as well as heat transfer due to the combined process of 
conduction and radiation. The program provides the general solution as well 
as the limiting (optically thin and large path length) solutions. The 
program is written for the five-band H 2 O gas, but can be used for any gas by 
employing appropriate thermophysical and radiative properties. 


D-1 
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1 PROGRAM TRAN (0UTPUT=65,TAPE6, TAPE 10, TAPE 15) 


THIS PROGRAM CALCULATES THE TRANSIENT RADIATIVE 
INTERACTION IN H20 BETWEEN TWO PARALLEL PLATES 
FOR CONSTANT WALL TEMPRATURE. 


T NONDIMENSIONAL TIME 

p PRESSURE (ATM) 

TW TEMPERATURE OF LOWER PLATE (K) 

L DISTANCE BETWEEN TWO PLATES (CM) 

Z NONDIMENSIONAL Y COORDINATE 

TM CHARACTERS ITIC TIME (SEC) 

IR=1 ONLY RADIATION 

IR=2 RADIATION & CONDUCTION 

IR=3 ONLY CONDUCTION 

EXTERNAL F311,F312,F313,F314,F315,F321,F322,F323,F324,F325,F411 
&,F412,F413,F414,F415,F421,F422,F423,F424,F425,FQ11,FQ12,FQ13,FQ14 
& ,FQ15,FQ21,FQ22,FQ23,FQ24,FQ25 
DIMENSION U(20) ,EPS(3) ,PRES( 10) ,TEMP(8) ,TT(20) 

REAL KFB,L,M1,M2,M3,M4,M5,R 

COMMON F1,F2,F3,F4,F5,B1,B2,B3,B4,B5,THETA,Z 

DATA U/0. 1,0. 2, 0.5, 1.0, 2., 5., 10., 20., 50., 100., 

& 200., 500., 1000./ 

DATA PRES/0.01,0.1,1.0,10.0/ 

DATA TEMP/300 .0,500.0, 1000 . 0 , 2000 . 0 , 3000 .0/ 

DATA TT/0.01,0.02,0.05,0. 1,0.5, 1.0/ 

DO 11 IR=1,3 


IF(IR.EQ.l) WRITE(6,111) 
IF(IR.EQ.l) WRITE(10,114) 
IF(IR.EQ.l) WRITE(15,115) 
IF(IR.EQ.2) WRITE(6,112) 
IF(IR.EQ.2) WRITE(10,113) 
IF(IR.EQ.2) WRITE(15,116) 
IF(IR.EQ.3) WRITE(6,561) 
IF(IR.EQ.3) WRITE(10,564) 
IF(IR.EQ.3) WRITE(15,565) 
TM=0. 00001 


WRITE(6,886) TM 
WRITE(10,886)TM 
WRITE(15,886)TM 


TIME MARCHING 

00 22 ITT=1,5 

T=TT(ITT) 

WRITE(6,888)T 

WRITE(10,888)T 

WRITE(15,888) T 

WRITE(6,666) 

WRITE(10,666) 

WRITE(15,666) 

00 33 IT=1,4 
TW=TEMP(IT) 



o o o o o o o o o o o o r> 
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T2*2*TW 

OZ=0.1 

Z=0.1 

DO 44 IZ=1,5 
DO 55 KK=2,4 
P=PRES(KK) 

DO 66 1=1,13 
L=U(I) 

IF(IR.EQ.3) GOTO 360 

CALCULATION OF PLANCK FUNCTION AND ITS DERIVATIVE 

WNB BAND CENTER (1/CM) 

HCK CONSTANT (K CM) 

CCC C1*C2 (ERG-K-CM**3/SEC) 

PFDBI— PLANCK FUNCTION DERIVATIVE FOR I BAND 

CESS=TW**2.0 

STU=TW**0.5 

SNT=TW/273.0 

HCK=1. 439257246 

BT3=5 . 668* 10** ( -5 . ) *TW**3 . 

WN1=500.0 

WN2=1600.0 

WN3=3750.0 

WN4=5350.0 

WN5=7250.0 

C2C1=HCK*WN1 

C2C2=HCK*WN2 

C2C3=HCK*WN3 

C2C4=HCK*WN4 

C2C5=HCK*WN5 

CCC=0. 000053847734 

SPECTROSCOPIC PROPERTIES OF H20 
WNBI BAND CENTER (1/CM) 

HERE WE HAVE CONSIDERED ONLY FIVE BANDS(20,6.3,2.7,1.87, &1.38 MICRONS) 

CCC1=CCC*(WN1**4) 

CCC2=CCC*(WN2**4) 

CCC3=CCC*(WN3**4) 

CCC4=CCC*(WN4**4) 

CCC5=CCC*(WN5**4) 

TB1=C2C1/TW 

TB2=C2C2/TW 

TB3=C2C3/TW 

TB4=C2C4/TW 

TB5=C2C5/TW 

TEB1=EXP(TB1) 

TEB2=EXP(TB2) 

TEB3=EXP(TB3) 

TEB4=EXP(TB4) 

TEB5=EXP(TB5) 

SNTB1=CESS*((TEB1-1.0)**2.0) 

SNTB2=CESS*((TEB2-1.0)**2.0) 



o o o o o o 
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SNTB3=CESS*((TEB3-1.0)**2.0) 
SNTB4=CESS*((TEB4-1.0)**2.0) 
SNTB5=CESS*((TEB5-1.0)**2.0) 

PFDB1=(CCC1*TEB1)/SNTB1 
PFDB2=(CCC2*TEB2)/SNTB2 
PFDB3=(CCC3*TEB3)/SNTB3 
PFDB4=(CCC4*TEB4)/SNTB4 
PFDB5=(CCC5*TEB5)/SNTB5 

BAND MODEL CORRELATIONS (TIEN & LOWDER WIDE BAND MODEL) 

AZI AOI/(l/CM) 

CZSI COI**2 (1/(ATM-CM)) 

BSI B**2 (NON DIMENSIONAL) 

OMGI WAVE NUMBER (1/CM) 

AKl=(TW/300.0)**0.5 
AK2=(300.0/TW)**1.5 
AZ1=49.4*AK1 
AZ2*90.1*AK1 
AZ3=112.6*AK1 
AZ4=79.7*AK1 
AZ5*79.7*AK1 
0MG1=3652.0 
0MG2=1595.0 
0MG3=3756.0 
TX=-(HCK/TW) 

TX1=TX*0MG1 
TX2=TX*0MG2 
TX3=TX*0MG3 
ETX1=EXP(TX1) 

ETX2=EXP(TX2) 

ETX3=EXP(TX3) 

C1=1.0-ETX1 
C2=1.0-ETX2 
C3=1.0-ETX3 
BRKT=TX*(OMGl+OMG3) 

PHI1Z1=(1.0-EXP(BRKT))/(C1*C2*C3) 
DJ=17.6/(TW/100.0)**0.5 
PHI7=EXP(-DJ) 

TS*300.0/TW 
CZS1*771*AK2*PHI7 
CZS2=3.35*AK2 
CZS3=1.52*AK2 
CZS4=0.276*AK2*PHI1Z1 
CZS5=0.23*AK2*PHI1Z1 

SI INTEGRATED BAND INTENSITY (1/(ATM CM**2)) 

S1=CZS1*AZ1 
S2=CZS2*AZ2 
S3=CZS3*AZ3 
S4=CZS4*AZ4 
S5=CZS5*AZ5 
C 

C PL PRESSURE PATH LENGTH (ATM-CM) 


/ 


ooooo ooo ooo ooo 
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c 

PL=P*L 
UZ1=CZS1*PL 
UZ2=CZS2*PL 
UZ3=CZS3*PL 
UZ4=CZS4*PL 
UZ5=CZS5*PL 
B1=1.5*UZ1 
B2*1.5*UZ2 
B3=1.5*UZ3 
B4=1.5*UZ4 
B5=1.5*UZ5 
BS1=0.073/AK1 
BS2*0.13/AK1 
BS3=0.145/AK1 
BS4=0.118/AK1 
BS5=0.201/AK1 

PEI EFFECTIVE PRESSURE FOR EACH BAND (NON DIMENSIONAL) 

PE1=(5.0*P) 

PE2=(5.0*P) 

PE3=(5.0*P) 

PE4=(5.0*P) 

PE5=(5.0*P) 

beta I LINE STRUCTURE PARAMETER 

BETA1=BS1*PE1 
BETA2=BS2*PE2 
BETA3=BS3*PE3 
BETA4=BS4*PE4 
BETA5=BS5*PE5 

CORRELATION FOR EACH BAND 

F1=2.94*(1.0-EXP(-(2.60*BETA1))) 

F2=2 . 94* ( 1 . 0-EXP ( - ( 2 . 60*BETA2 ) ) ) 

F3=2 . 94* ( 1 . 0-EXP ( - ( 2 . 60*BETA3 ) ) ) 

F4=2 . 94* ( 1 . 0-EXP ( - { 2 . 60*BETA4) ) ) 
F5=2.94*(1.0-EXP(-{2.60*BETA5))) 

NUMERICAL INTEGRATION 

FOR DETAILS OF INTEGRATION REFER MATHEMATICAL LIBRARY 
AT NASA LANGLEY (DOCUMENT N-3) 

EPS(1)=1.0E-12 
EPS(2)=1.0E-12 
X=0.0 

CALL CADRE(X,Z,F311,EPS,ITEXT,G311,IERR) 

CALL CADRE(X,Z,F312,EPS,ITEXT,G312,IERR) 

CALL CADRE(X,Z,F313,EPS,ITEXT,G313,IERR) 

CALL CADRE(X,Z,F314,EPS,ITEXT,G314,IERR) 

CALL CADRE(X,Z,F315,EPS,ITEXT,G315,IERR) 

CALL CADRE(Z,1.,F321,EPS,ITEXT,G321,IERR) 
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CALL CADRE(Z,1.,F322,EPS,ITEXT,G322,IERR) 
call CA0RE(Z,1.,F323,EPS,ITEXT,G323,IERR) 

CALL CA0RE(Z,1.,F324,EPS,ITEXT,G324,IERR) 

CALL CADRE(Z,1.,F325,EPS,ITEXT,G325,IERR) 

CALL CAORE(X, Z,F411, EPS, ITEXT, 6411, lERR) 

CALL CADRE(X, Z,F412, EPS, ITEXT, G412,IERR) 

CALL CADRE(X, Z,F413, EPS, ITEXT, G413,IERR) 

CALL CAORE(X, Z,F414, EPS, ITEXT, G414,IERR) 

CALL CAORE(X, Z,F415, EPS, ITEXT, G415,IERR) 

CALL CADRE(Z, 1. 0,F421, EPS, ITEXT, G421,IERR) 

CALL CADRE ( Z , 1 . 0 , F422 , EPS , ITEXT , G422 , I ERR ) 

CALL CADRE(Z, 1. 0,F423, EPS, ITEXT, G423,IERR) 

CALL CAORE(Z, 1. 0,F424, EPS, ITEXT, G424,IERR) 

CALL CAORE(Z, 1.0,F425, EPS, ITEXT, G425,IERR) 

H1=AZ1*PFD81 
H2=AZ2*PFDB2 
H3=AZ3*PFDB3 
H4=AZ4*PFDB4 
H5=AZ5*PFDB5 
360 CONTINUE 

RHO DENSITY 

CP SPECIFIC HEAT AT CONSTANT PRESSURE. 

KFB THERMAL CONDUCTIVITY 

RH0=0. 0128677*18. /TW 
TC=TW/100. 

CP= ( 143 . 05- 183 . 54*TC**0 . 25+82 . 7 5 1*TC**0 . 5- 3 . 6989*TC ) / 18 . 
CON=TM/(L*RHO*CP) 

KFB=(1511.620876)*(SNT**1.48) 

R=(KFB*TM)/(L*L*RHO*CP) 

Z2=Z*Z 
CZ=1./(Z-Z2) 

IF(IR.EQ.3) GOTO 362 
M1=C0N*H1 
M2=C0N*H2 
M3=CON*H3 
M4=C0N*H4 
M5=C0N*H5 
AN1=M1*UZ1 
AN2=M2*UZ2 
AN3=M3*UZ3 
AN4=M4*UZ4 
AN5=M5*UZ5 

SUM31=AN1*1.5*(G311-G321) 

SUM32=AN2*1 . 5*(G312-G322) 

SUM33=AN3*1.5*(G313-G323) 

SUM34=AN4*1 . 5*(G314-G324) 

SUM35=AN5*1.5*(G315-G325) 
63=CZ*(SUM31+SUM32+SUM33+SUM34+SUM35) 

IF(IR.EQ.l) R=0.0 
AJ3=G3+CZ*2.*R 
SUM41=AN1*(G411-G421) 

SUM42=AN2*(G412-G422) 

SUM43=AN3*(G413-6423) 
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SUM44=AN4*(G414-G424) 

SUM45*AN5*(G415-G425) 

G4=-3*CZ* ( SUM41+SUM42+SUM43+SUM44+SUM45 ) 

AJ4=G4+2.*CZ*R 

GZ=(1-Z2)*CZ 

EAJ3T=EXP(-(AJ3*T)) 

GT=(GZ-(AJ4/AJ3) )*EAJ3T+AJ4/AJ3 
362 CONTINUE 

IF(IR.EQ.3) GT=(EXP(2*R*T/(Z2-Z))+Z)/Z 
THETA=Z2+GT*(Z-Z2) 

IF(IR.EQ.3) DFLUX=-KFB*(T2-TW)*(2*Z+GT*(1-2*Z))/L 

TW4=TW**4 

T24=T2**4 

I F ( I R . EQ . 3 ) FLUX=OFLUX/ ( 5 . 668* 10** ( -5 . ) * ( TW4-T24 ) ) 

IF(IR.EQ.3) GOTO 364 

CALL CADRE(0.0,Z,FQ11,EPS,ITEXT,GQ11.IERR) 

CALL CADRE(0.0,Z,FQ12,EPS,ITEXT,GQ12,IERR) 

CALL CADRE(0.0,Z,FQ13,EPS,ITEXT,GQ13JERR) 

CALL CADRE(0.0,Z,FQ14,EPS,ITEXT,GQ14,IERR) 

CALL CADRE(0.0,Z,FQ15,EPS,ITEXT,GQ15,IERR) 

CALL CADRE(Z,1.0,FQ21,EPS,ITEXT,GQ21,IERR) 

CALL CA0RE(Z,1.0,FQ22,EPS,ITEXT,GQ22,IERR) 

CALL CA0RE(Z,l.O,FQ23,EPS,ITEXT,GQ23,IERR) 

CALL CADRE(Z,1.0,FQ24,EPS,ITEXT,GQ24,IERR) 

CALL CA0RE(Z,1.0,FQ25,EPS,ITEXT,GQ25,IERR) 

FLUX=1-3.*(UZ1*H1*(GQ11+GQ21)+UZ2*H2*(GQ12+GQ22)+UZ3*H3*(GQ13+GQ23 

&)+UZ4*H4*(GQ14+GQ24)+UZ5*H5*(GQ15+GQ25))/(8*BT3) 

DFLUX=5,668*10**(-5.)*FLUX*(TW4-T24) 

364 CONTINUE 

WRITE(6,222) TW,P,L, THETA, Z, FLUX, DFLUX 
IF(IR.EQ.3) GOTO 366 

LARGE PATH LENGTH 

DS=ABS((Z-1)/Z) 

G 1 5= ( 2- ( 2*Z- 1 ) *ALOG ( OS ) ) *CZ* ( M 1+M2+M3+M4+M5 ) 

G16= ( 1-Z*AL0G (OS ) ) *CZ*(M1+M2+M3+M4+M5 ) *2 . 

AJ15=G15+2.*R*CZ 
AJ16=G16+2.*R*CZ 
EAJ15=EXP(-AJ15*T) 

G4T=(GZ-(AJ16/AJ15))*EAJ15+AJ16/AJ15 
Z2=Z*Z 

THEL=Z2+G4T*(Z-Z2) 

GQ3=THEL*(AL0G(Z)+13. 81551055) 

GQ4=(1-THEL)*(AL0G(1-Z)+13. 81551055) 
FLUXL=1-(H1+H2+H3+H4+H5)*(GQ3+GQ4)/(4*BT3) 

WRITE(10,222) TW,P,L,THEL,Z ,FLUXL 

OPTICALLY THIN — H20 

AJ25=3*(AN1+AN2+AN3+AN4+AN5)+2*R/(Z-Z2) 
AJ26=(-3*(ANl+AN2+AN3+AN4+AN5)*(Z2-0.5)+2*R)/(Z-Z2) 
EAJ25=EXP(-(AJ25*T)) 

G6T= (GZ- ( AJ26/AJ25) ) *EAJ25+AJ26/AJ25 
THET=Z2+G6T*(Z-Z2) 
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SPF0=S1*PF0B1+S2*PFDB2+S3*PFDB3+S4*PFDB4+S5*PFD85 

FLUX0=l-(3./(8*BT3))*P*L*SPFD*((l-Z)+(2*Z-l)*THET) 

WRITE( 15,222) TW,P,L,THET,Z,FLUXO 
366 CONTINUE 

66 CONTINUE 

55 CONTINUE 

Z=Z+OZ 

44 CONTINUE 

33 CONTINUE 

22 CONTINUE 

II CONTINUE 

III FORMAT(5X,*GENERAL — H20 ONLY RADIATION*) 

112 FORMAT(5X,*GENERAL — H20 RADIATION AND CONDUCTION*) 

113 FORMAT (5X,*LARGE PATH LENGTH — H20 RADIATION AND CONDUCTION*) 

114 FORMAT (5X,*LARGE PATH LENGTH — H20 ONLY RADIATION*) 

115 FORMAT(5X,* OPTICALLY THIN — H20- — ONLY RADIATION*) 

116 FORMAT(5X,* OPTICALLY THIN — H20 RADIATION AND CONDUCTION*) 

561 FORMAT(5X,*GENERAL — H20 ONLY CONDUCTION*) 

564 FORMAT (5X,*LARGE PATH LENGTH — H20 ONLY CONDUCTION*) 

565 FORMAT(5X,* OPTICALLY THIN — H20 ONLY CONDUCTION*) 

886 F0RMAT(10X,*CHARACTERSITIC TIME (TM)=*,F12.9) 

888 F0RMAT(10X,*TIME=*,F8.3) 

666 F0RMAT(4X,*TW*,T18,*P*,T28,*L*,T38,*THETA*,T48,*Z*,T56,*FLUX* 

&,T65,*DFLUX*) 

222 F0RMAT(6(1X, F9.4) ,1X,E12.6) 

STOP 

END 

FUNCTION F311(ZB) 

COMMON F1,F2,F3,F4,F5,B1,B2,B3,B4,B5,THETA,Z 
XX=B1*(Z-ZB) 

DENl=(Fl*((XX**2.0)+(2.0*XX)+2.0)+XX)*(XX+2.0*Fl) 

AUD1=(F1*((XX**2.0)+(4.0*XX*F1)+(4.0*F1)))/DEN1 

F311=(1-2*ZB)*AUD1 

RETURN 

END 

FUNCTION F312(ZB) 

COMMON F1,F2,F3,F4,F5,B1,B2,83,B4,B5,THETA,Z 
XX=B2*(Z-ZB) 

OEN2=(F2*((XX**2.0)+{2.0*XX)+2.0)+XX)*(XX+2.0*F2) 

AUD2*(F2*{(XX**2.0)+(4.0*XX*F2)+(4.0*F2)))/0EN2 

F312=(1-2*ZB)*AUD2 

RETURN 

END 

FUNCTION F313(ZB) 

COMMON F1,F2,F3,F4,F5,B1,B2,B3,B4,85,THETA,Z 
XX=B3*(Z-ZB) 

DEN3=(F3*((XX**2.0)+(2.0*XX)+2.0)+XX)*(XX+2.0*F3) 

AUD3=(F3*((XX**2.0)+(4.0*XX*F3)+(4.0*F3)))/DEN3 

F313=(1-2*ZB)*AUD3 

RETURN 

END 

FUNCTION F314(ZB) 

COMMON FI, F2.F3,F4,F5.B1,B2,B3,B4,B5, THETA, Z 
XX=B4*(Z-ZB) 

OEN4=(F4*((XX**2.0)+(2.0*XX)+2.0)+XX)*(XX+2.0*F4) 



AUD4=(F4*((XX**2.0)+(4.0*XX*F4)+(4.0*F4)))/DEN4 

F314=(1-2*ZB)*AUD4 

RETURN 

END 

FUNCTION F315(ZB) 

COMMON F1,F2,F3,F4,F5,B1,B2,B3,B4,B5,THETA,Z 
XX=B5*(Z-ZB) 

DEN5=(F5*((XX**2.0)+(2.0*XX)+2.0)+XX)*(XX+2.0*F5) 

AUD5=(F5*((XX**2.0)+(4.0*XX*F5)+(4.0*F5)))/DEN5 

F315=(1-2*ZB)*AUD5 

RETURN 

END 

FUNCTION F321(ZB) 

COMMON FI, F2,F3,F4,F5,B1,B2,B3,B4,B5, THETA, Z 
XX=B1*(ZB-Z) 

OEN1=(F1*((XX**2.0)+(2.0*XX)+2.0)+XX)*(XX+2.0*F1) 

AUD1*(F1*((XX**2.0)+(4.0*XX*F1)+(4.0*F1)))/DEN1 

F321=(1-2*ZB)*AUD1 

RETURN 

END 

FUNCTION F322(ZB) 

COMMON F1,F2,F3,F4,F5,B1,B2,B3,84,B5,THETA,Z 
XX=B2*(ZB-Z) 

DEN2=(F2*((XX**2.0)+(2.0*XX)+2.0)+XX)*(XX+2.0*F2) 

AUD2»(F2*((XX**2.0)+(4.0*XX*F2)+(4.0*F2)))/DEN2 

F322=(1-2*ZB)*AU02 

RETURN 

END 

FUNCTION F323(ZB) 

COMMON FI, F2,F3,F4,F5,B1,B2,B3,B4,B5, THETA, Z 
XX=B3*(ZB-Z) 

DEN3=(F3*((XX**2.0)+{2.0*XX)+2.0)+XX)*(XX+2.0*F3) 

AUD3=(F3*((XX**2.0)+(4.0*XX*F3)+(4.0*F3)))/DEN3 

F323=(1-2*ZB)*AUD3 

RETURN 

END 

FUNCTION F324(ZB) 

COMMON F1,F2,F3,F4,F5,B1,B2,B3,B4,B5,THETA,Z 
XX=B4*(ZB-Z) 

OEN4=(F4*((XX**2.0)+(2.0*XX)+2.0)+XX)*(XX+2.0*F4) 

AUD4=(F4*((XX**2.0)+(4.0*XX*F4)+(4.0*F4)))/DEN4 

F324=(1-2*ZB)*AUD4 

RETURN 

END 

FUNCTION F325(ZB) 

COMMON FI, F2,F3,F4,F5,B1,B2,B3,B4,B5, THETA, Z 
XX=B5*(ZB-Z) 

DEN5=(F5*((XX**2.0)+(2.0*XX)+2.0)+XX)*(XX+2.0*F5) 

AUD5=(F5*((XX**2.0)+(4.0*XX*F5)+(4.0*F5)))/DEN5 

F325=(1-2*ZB)*AUD5 

RETURN 

END 

FUNCTION F411(ZB) 

COMMON FI, F2,F3,F4,F5,B1,B2,B3,B4,B5, THETA, Z 
XX=81*(Z-ZB) 



DENl=(Fl*((XX**2.0)+{2.0*XX)+2.0)+XX)*(XX+2.0*Fl) 

AUD1=(F1*((XX**2.0)+(4.0*XX*F1)+(4.0*F1)))/DEN1 

F411=ZB*AU01 

RETURN 

END 

FUNCTION F412(ZB) 

COMMON FI, F2,F3,F4,F5,B1,B2,B3,B4,B5, THETA, Z 
XX=B2*(Z-ZB) 

DEN2=(F2*((XX**2.0)+(2.0*XX)+2.0)+XX)*(XX+2.0*F2) 

AU02=(F2*((XX**2.0)+(4.0*XX*F2)+(4.0*F2)))/OEN2 

F412=ZB*AUD2 

RETURN 

END 

FUNCTION F413(ZB) 

COMMON FI, F2,F3,F4,F5,B1,B2,B3,B4,B5, THETA, Z 
XX=B3*(Z-ZB) 

DEN3=(F3*((XX**2.0)+(2.0*XX)+2.0)+XX)*(XX+2.0*F3) 

AU03=(F3*((XX**2.0)+(4.0*XX*F3)+{4.0*F3)))/DEN3 

F413*ZB*AU03 

RETURN 

END 

FUNCTION F414(ZB) 

COMMON FI, F2,F3,F4,F5,B1,B2,B3,B4,B5, THETA, Z 
XX=B4*(Z-ZB) 

DEN4=(F4*((XX**2.0)+(2.0*XX)+2.0)+XX)*(XX+2.0*F4) 

AU04=(F4*( (XX**2.0)+(4.0*XX*F4)+(4.0*F4) ) )/DEN4 

F414=ZB*AUD4 

RETURN 

END 

FUNCTION F415(ZB) 

COMMON FI, F2,F3,F4,F5,B1,B2,B3,B4,B5, THETA, Z 
XX=B5*(Z-ZB) 

DEN5=(F5*((XX**2.0)+(2.0*XX)+2.0)+XX)*(XX+2.0*F5) 

AUD5=(F5*((XX**2.0)+(4.0*XX*F5)+(4.0*F5)))/DEN5 

F415=ZB*AUD5 

RETURN 

END 

FUNCTION F421(ZB) 

COMMON FI, F2,F3,F4,F5,B1,B2,B3,B4,B5, THETA, Z 
XX=B1*(ZB-Z) 

DENl=(Fl*((XX**2.0)+(2.0*XX)+2.0)+XX)*(XX+2.0*Fl) 

AUD1=(F1*((XX**2.0)+(4.0*XX*F1)+(4.0*F1)))/DEN1 

F421*ZB*AUD1 

RETURN 

END 

FUNCTION F422(ZB) 

COMMON FI, F2,F3,F4,F5,B1,B2,B3,B4,B5, THETA, Z 
XX=B2*(ZB-Z) 

DEN2=(F2*((XX**2.0)+(2.0*XX)+2.0)+XX)*(XX+2.0*F2) 

AUD2=(F2*({XX**2.0)+(4.0*XX*F2)+(4.0*F2)))/0EN2 

F422=ZB*AUD2 

RETURN 

END 

FUNCTION F423(ZB) 

COMMON FI, F2,F3,F4,F5,B1,B2,B3,B4,B5, THETA, Z 



XX=B3*(ZB-Z) 

DEN3=(F3*((XX**2.0)+(2.0*XX)+2.0)+XX)*(XX+2.0*F3) 

AUD3=(F3*((XX**2.0)+{4.0*XX*F3)+(4.0*F3)))/DEN3 

F423=ZB*AU03 

RETURN 

END 

FUNCTION F424(ZB) 

COMMON FI, F2,F3,F4,F5,B1,B2,B3,B4,B5, THETA, Z 
XX=B4*(ZB-Z) 

0EN4»(F4*((XX**2.0)+{2.0*XX)+2.0)+XX)*(XX+2.0*F4) 

AUD4= ( F4* ( ( XX**2 . 0 )+{ 4 . 0*XX*F4)+{ 4 . 0*F4 ) ) ) /DEN4 

F424=ZB*AU04 

RETURN 

END 

FUNCTION F425(ZB) 

COMMON FI, F2,F3,F4,F5,B1,B2,B3,B4,B5, THETA, Z 
XX=B5*(ZB-Z) 

DEN5*(F5*((XX**2.0)+(2.0*XX)+2.0)+XX)*(XX+2.0*F5) 

AUD5=(F5*((XX**2.0)+(4.0*XX*F5)+(4.0*F5)))/DEN5 

F425=ZB*AU05 

RETURN 

END 

FUNCTION FQll(ZB) 

COMMON F1,F2,F3,F4,F5,B1,B2,B3,84,B5,THETA,Z 
XX=B1*(Z-ZB) 

DENl=(Fl*((XX**2.0)+(2.0*XX)+2.0)+XX)*(XX+2.0*Fl) 

AUD1=(F1*((XX**2.0)+(4.0*XX*F1)+{4.0*F1)))/DEN1 

FQ11=(THETA)*AUD1 

RETURN 

END 

FUNCTION FQ12(ZB) 

COMMON FI, F2,F3,F4,F5,B1,B2,B3,B4,B5, THETA, Z 
XX=B2*(Z-ZB) 

DEN2=(F2*((XX**2.0)+(2.0*XX)+2.0)+XX)*(XX+2.0*F2) 

AUD2=(F2*((XX**2.0)+(4.0*XX*F2)+(4.0*F2)))/DEN2 

FQ12=(THETA)*AUD2 

RETURN 

END 

FUNCTION FQ13(ZB) 

COMMON FI, F2,F3,F4,F5,B1,B2,B3,B4,B5, THETA, Z 
XX=B3*(Z-ZB) 

DEN3=(F3*((XX**2.0)+(2.0*XX)+2.0)+XX)*(XX+2.0*F3) 

AUD3=(F3*((XX**2.0)+(4.0*XX*F3)+(4.0*F3)))/DEN3 

FQ13=(THETA)*AUD3 

RETURN 

END 

FUNCTION FQ14(ZB) 

COMMON FI, F2,F3,F4,F5,B1,B2,B3,B4,B5, THETA, Z 
XX=B4*(Z-ZB) 

D£N4=(F4*((XX**2.0)+(2.0*XX)+2.0)+XX)*(XX+2.0*F4) 

AUD4=(F4*((XX**2.0)+(4.0*XX*F4)+(4.0*F4)))/DEN4 

FQ14=(THETA)*AUD4 

RETURN 

END 

FUNCTION FQ15(ZB) 



COMMON F1,F2,F3,F4,F5,B1,B2,B3,B4,85,THETA,Z 
XX=B5*(Z-ZB) 

DEN5*(F5*((XX**2.0)+(2.0*XX)+2.0)+XX)*(XX+2.0*F5) 

AUD5=(F5*((XX**2.0)+{4.0*XX*F5)+(4.0*F5)))/DEN5 

FQ15=(THETA)*AUD5 

RETURN 

END 

FUNCTION FQ21(ZB) 

COMMON FI, F2,F3,F4,F5,B1,B2,B3,B4,B5, THETA, Z 
XX=B1*(ZB-Z) 

DEN1=(F1*((XX**2.0)+(2.0*XX)+2.0)+XX)*(XX+2.0*F1) 

AUD1=(F1*( (XX**2.0)+(4.0*XX*F1)+(4.0*F1) ) )/DENl 

FQ21=( 1-THETA) *AUD1 

RETURN 

END 

FUNCTION FQ22(ZB) 

COMMON FI, F2,F3,F4,F5,B1,B2,B3,B4,B5, THETA, Z 
XX=B2*(ZB-Z) 

OEN2=(F2*((XX**2.0)+(2.0*XX)+2.0)+XX)*(XX+2.0*F2) 

AUD2=(F2*((XX**2.0)+(4.0*XX*F2)+(4.0*F2)))/DEN2 

FQ22=( 1-THETA) *AUD2 

RETURN 

END 

FUNCTION FQ23(ZB) 

COMMON FI, F2,F3,F4,F5,B1,B2,B3,B4,B5, THETA, Z 
XX=B3*(ZB-Z) 

DEN3=(F3*((XX**2.0)+(2.0*XX)+2.0)+XX)*(XX+2.0*F3) 

AUD3=(F3*((XX**2.O)+(4.O*XX*F3)+(4.0*F3)))/OEN3 

FQ23=(1-THETA)*AUD3 

RETURN 

END 

FUNCTION FQ24(ZB) 

COMMON FI, F2,F3.F4,F5,B1,B2,B3,B4,B5, THETA, Z 
XX=B4*(ZB-Z) 

DEN4=(F4*((XX**2.0)+(2.0*XX)+2.0)+XX)*(XX+2.0*F4) 

AUD4=(F4*((XX**2.0)+{4.0*XX*F4)+(4.0*F4)))/DEN4 

FQ24=(1-THETA)*AUD4 

RETURN 

END 

FUNCTION FQ25(ZB) 

COMMON FI, F2,F3,F4,F5,B1,B2,B3,B4,B5, THETA, Z 
XX=B5*(ZB-Z) 

DEN5=(F5*((XX**2.0)+(2.0*XX)+2.0)+XX)*(XX+2.0*F5) 

AUD5*(F5*((XX**2.0)+(4.0*XX*F5)+(4.0*F5)))/DEN5 

FQ25=(1-THETA)*AU05 

RETURN 

END 


